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Abstract 



We report on a study of a classical, finite system of confined particles 
in two dimensions with a two-body repulsive interaction. We first develop a 
simple analytical method to obtain equilibrium configurations and energies 
for few particles. When the confinement is harmonic, we prove that the first 
transition from a single shell occurs when the number of particles changes 
from five to six. The shell structure in the case of an arbitrary number of 
particles is shown to be independent of the strength of the interaction but 
dependent only on its functional form. It is also independent of the magnetic 
field strength when included. We further study the effect of the functional 
form of the confinement potential on the shell structure. Finally, we report 
some interesting results when a three-body interaction is included, albeit in 
a particular model. 
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I. INTRODUCTION 



Two dimensional clusters or "artificial atoms" have attracted considerable attention in 
recent years Jl||. A cluster of a finite number of charged particles confined by an external 
potential may be regarded as an "artificial atom" . There are several examples of such systems 
from mesoscopic systems to the astroplasma system. In observed systems such as electrons 
on a liquid helium film f|j , drops of colloidal suspensions [[| and confined dusty particles [£| , 
the dynamics is essentially classical. On the other hand, in mesoscopic systems like quantum 
dots, quantum effects may not be negligible || . The actual Hamiltonians which capture the 
full dynamics of these systems are rather complicated. However, several studies assuming 
model systems have been carried out. In particular, there have been studies on the ordering 
and transitions of charged particles in two dimensions pHTU]. 

In a recent Monte Carlo study, Bedanov and Peeters (see also Bolten and Rossler 
pi|) have analysed the classical ground state of a system of confined, charged particles 
interacting through the Coulomb interaction. By minimising the classical energy, they obtain 
numerically the shell structure in a cluster of A particles. They have systematically listed the 
shell structure in a "Mendeleev" table for A < 52, and for a few large clusters. Similar results 
are available also for logarithmic two-body interaction ||. An important fact that follows 
from this Mendeleev table is that the charged particles, confined in a parabolic potential, 
arrange themselves in concentric shells where each shell may be thought of as an annulus 
whose width is much smaller than the radius. For particular numbers, A = 6,19,37,... 
these annuli are almost precisely circles. These are called magic numbers. There are further 
systematics (approximate) in this table. Up to A = 5, the particles in the ground state 
configuration arrange themselves on a circle. For A = 6, there are five particles on the circle 
and one in the centre (circle-dot) in the ground state. With the addition of more particles, 
a second shell starts forming until the arrangement has 5 particles in the inner shell and 10 
particles on the outer shell, ie A = 15. Addition of one more particle (A = 16), creates a 
configuration similar to A = 15 with the extra particle at the centre. It is interesting to 
note that the multi-shell structure was experimentally observed by Christian Myer | I 



more 



than a hundred years ago. He observed these geometric transitions in a system involving 
magnetic needles floating on water, confined by a bar magnet held above the surface. The 
result, as quoted by J. J. Thompson, notes that five magnetic needles always remained on the 
circumference of a circle whereas the sixth one, when added, always drifted to the middle. 
Notice that the interaction between magnetic needles is different from Coulomb interaction. 
Nevertheless the observed phenomenon bears resemblance to the numerical results. 

At least up to A = 50, the geometric transition from 5 to 6, that is from a circle to 
a circle-dot, forms a template for the innermost shell. To a very good approximation, the 
number of shells for any A, in the Mendeleev table ||, may be deduced as follows: Write a 
given A as a sum of non-repeating multiples of five and a remainder. The number of terms in 
the sum is the number of shells. For example, A = 48 may be written as A = 5 + 10 + 15 + 13. 
This, in our way of looking, corresponds to having four shells. (It is not expected that this 
will hold for large A since the system goes into a hexagonal arrangement in the bulk.) The 
main theme of the paper is this geometric transition from circle to circle-dot configuration 
when the number of particles changes from 5 to 6. We also derive some general results 
concerning the shell- structure for arbitrary A. 
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In this paper we first consider a cluster of N particles in which the repulsive two-body 
potential is a power-law. The classical system we are interested in, to begin with, consists of 
N particles of equal mass m, confined in an oscillator potential, in a uniform magnetic field 
and interacting via a two-body potential. The Hamiltonian of such a system of particles is 
given by 
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where and denote the position and momentum vectors of the particle with index i. The 
vector potential, for a uniform magnetic field in the symmetric gauge, is given by 

{a%)x = -^lVi ; (ai) y = uj L Xi, (2) 

where uo^ is the Larmor frequency and = fi — fj. The power v (positive) is kept arbitrary. 
The Coulomb Hamiltonian is recovered by choosing v = 1/2. We also consider the case 
when the two-body interaction is of the form — j3 Ylijtyi) log^/p 2 ) which is repulsive for 
rfj < p 2 (see Appendix A). This may be closer to the real situation as also the v — 1/2 



(Coulomb) case for electrons in a quantum dot ||12[| . Recently, the v — 1 (inverse square 
interaction) case has also been analysed in detail |13|] because of its relevance to quantum 
dot systems. In addition, this case nicely lends itself to analytical manipulations. 

We devise a simple analytical method to obtain the classical ground state energy in 
two steps. First, we extremise the energy for a fixed total angular momentum J (which 
is conserved), and then minimise this energy with respect to J. This has two advantages. 
It reflects the quantum degeneracy of the lowest Landau level for electrons in a uniform 
magnetic field in the absence of any interaction, at the classical level itself and secondly it 
allows one to do the second step of the minimisation over the quantised values of the angular 
momentum. Here, however, we restrict ourselves to a classical analysis only and derive some 
exact results analytically for the equilibrium configurations. 

In Section II, we show that: 

1. The configurations extremising the energy, modulo overall scale, are independent of 
the parameters of the Hamiltonian and also of the angular momentum J so long as 
the repulsive two-body potential falls off as a power-law with the relative distance and 
the confinement is harmonic. Only the overall length scale is sensitive to these details. 

2. Two special configurations, one in which all the N particles are on a circle (referred 
to as O) an d the one in which N — 1 are on a circle with one at the center (referred 
to as O) are always equilibrium configurations. We give exact analytical expressions 
for the corresponding energy for all N and v. The O has lower energy for N < 5 
while O has lower energy for N > 6. In fact the O is the ground state for all iV < 5 
while O is the ground state for N = 6,7, 8. (This is a well known result for Coulomb 
interaction [T^H] but is also true in the general case considered below for arbitrary v 



up to a maximum value which depends on the number of particles.) This geometric 
transition in the ground state is the first one to occur and is independent of the precise 
form of the repulsive interaction. 
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3. While it is known numerically that for N > 9 and for Coulomb potential §JT^] the 
minimum energy configurations exhibit approximate multi-shell structure, the special 
configurations provide an upper bound on the minimum energy for a whole class of 
interactions that we consider here. The results presented in this section are an elabo- 



ration of our previous work [TJj] 



In Section III, we consider the effect of confinement on the geometric transition which 
occurs for N from 5 to 6. The one-body confinement in eq.([l]) is generalised to V con f = 
Y,iLi( r iV- The oscillator confinement is recovered for 7 = 1. Surprisingly, we find that 
there exists a critical value of v above which the first geometric transition always occurs 
between 5 and 6 for 7 > 1 and between 4 and 5 for 7 < 1. 

In Section IV, we consider the effect of three-body perturbations on this geometric tran- 
sition. We do so by using a particular model Hamiltonian whose exact quantum mechanical 
ground state energy and wave function are known in a limit to be defined later. While this 
again has interesting properties in its own right, it is included here primarily to study its 
effect on shell formation. Some of questions relating to three-body perturbations are being 
studied and will be published else where [|0 . 



We conclude in Section V with some comments and future prospects. Appendix A 
contains some general results concerning the logarithmic interaction potential. Details of 
our numerical simulations are discussed in Appendix B. These are used to check our results 
with ref. for the case of Coulomb potential and then extended to other forms of the 
potential. 



II. GEOMETRY OF CLUSTERS IN PARABOLIC CONFINEMENT 

In this section we consider the classical equilibrium configurations of the Hamiltonian 
given by eq. ([!]). The Hamiltonian can be written in terms of dimensionless units by intro- 



ducing a length scale / = y(h/(muj) which is the basic oscillator length. All distances are 
measured in terms of this basic length unit. Note that h is introduced only as a convenience 
so that the energy is measured in units of tkv and does not have any other significance as in 
the quantum case. The analysis presented in this paper is entirely classical. The momenta 
are measured in units oih/l. 

The new Hamiltonian in these scaled units, but keeping the same notation, may be 
written as 



TT N 



+ 9 E T^Tw ( 3 ) 



where ji = XiPi y — yiPi X , a — — and g = ^(l) 2l/ - Unless otherwise mentioned the summa- 
tions run from 1 to A hereafter. While the original coupling constant (3 was dimensional the 
new coupling constant g is dimensionless. Hereafter we assume all the energies are measured 
in units of fiui and do not write the units explicitly. 

To find the equilibrium configurations we carry out the variation in two steps. For the 
first step of the variation we introduce the function 

F = H + \(£ji-J), (4) 
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where ji are the single particle angular momenta and A is the Lagrange multiplier which en- 
forces the constraint J = Setting SF = 0, where the variation is done in the full phase 
space variables, gives the necessary equations to determine the equilibrium configuration in 
the phase space, 

Pix = (a + \)Vi, (5) 
p iy = -{a + X)xi, (6) 

(1 - A 2 - 2aA)3 = Agv £ -fy—. (7) 

This is the basic set of equations. Any solution to this set of equations describes an equi- 
librium configuration which could be a local minimum/ maximum or a saddle point. We 
also remark that while the Hamiltonian is written in a particular symmetric gauge, the 
variational equations given above are actually gauge invariant. The main advantage of in- 
troducing F instead of H for variation is that it allows us to keep the dependence of the 
variational equations on the magnetic field even at the classical level. If we vary H separately 
this advantage is lost. We may also include the case of logarithmic interaction by simply 
setting the power v = and by setting the prefactor to 4g instead of Agv in eq.(^). To 
avoid confusion, however, we briefly discuss the results for the logarithmic case separately 
in Appendix A. 

First we present a qualitative but general analysis of this basic set of equations without 
making any assumptions. To make the analysis simple, we introduce an auxiliary variable, 

* = (8) 

i 

The total angular momentum may now be written in terms of this auxiliary variable as 

J = ^2( x iPiy - ViVix) = + X)(f), (9) 

i 

where we have made use of eqs.(|5]||). It is convenient to express = Rsi where R is a 
common scale factor which may be taken to be the radius of the farthest particle, say the 
N th one, and Sj denotes the internal variables. Therefore 

JV-l 

<P = R 2 {Y / s* + l]=R 2 <P (10) 

i=l 

since sfy = 1 by choice. Using eq.flT]) and eliminating A-dependence using eq.(H), we have 

^ (s 2 - 



[l + a 2_^_^ J2 -i^-r. (11) 



Taking the scalar product with Sj and dividing both sides by s 2 0), we get 

(R 2 r +1 n | * _ = v i - (sj/ Si ) cos(%) 

*9" 1 R*r ^ l) {s 2 + s 2 -2s l s 3C o^6 tj )Y^ 1 ) 
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Note that the LHS is independent of the particle index i. Thus we have N — 1 independent 
of equations of the type 



E 



(sj/si) cos(%) 



Further by taking the cross product with Si and dividing by Sj, we get, 

0, (14) 



which provides a further set of N conditions on the internal coordinates Sj. Notice that 
these conditions are manifestly scale invariant. Together, eqs.(|TB|JT^) provide the 2N — 1 
necessary equations for determining the Sj and the angles 0j. These determining equa- 
tions are also completely independent of a, J and g. We have, therefore, the result that 
(si, S2, Sjv-i, #2, •••) #/v) are independent of the magnetic field (a), the total angular 
momentum J and the interaction strength g. These parameters, however, determine the 
overall scale R through eq.([L2|). In fact since s% = 1, the corresponding equation may be 
taken to be the determining equation for R in terms of the parameters of the Hamiltonian, 



kgv R 4 (p 2 

N-l 

W+l " 



v l- Sj cos{6 Nj ) 
3 tl (l + s]-2 Sj cos(6 Nj )y+i- 

Thus we have a very general result that the geometry, that is, the configuration modulo a 
common scale factor, or the shell structure of the equilibrium configuration is independent of 
the parameters of the Hamiltonian. The overall size of the system depends on the external 
magnetic field strength as well as on the total angular momentum J. (For further details 
on the influence of the magnetic field on the inter-shell rotation and diffusion, see ref. [[T5|]). 
The shell structure, however, depends on the nature of the repulsive interaction through 
the parameter v but not its strength. The above analysis is valid even if any one of the 
Si is zero, i.e., one particle being at the origin of the coordinate system, in which case we 
have two equations less (not more than one particle can be at the origin). Note that the 
statement above applies to all the extremum configurations and is not just restricted to the 
local minima . In Fig.l we have shown the typical shell structure for the ground state of 25 
particles as a function of J and v. It is easy to see that the shell structure, i.e., the number 
of particles in each shell and their relative orientation, is unchanged with J (modulo overall 
rotation). The shells get distorted with changes in v. Note that the system size has been 
normalised to the same value in Fig.l for convenience in display. 

The energy of the equilibrium configuration can be easily computed by noting that the 
auxiliary variable <fi defined in eq.(H) is related to the two-body potential energy by 



:i + ^-^)4> = 2gu £ * (16) 
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where the RHS is proportional to the potential energy due to interaction. Since the RHS 
and are positive definite we have the condition (1 + a 2 ) (ft 2 > J 2 . 
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Figure 1: Shell structure of the ground state for 25 particles as a function of the total 
angular momentum J and v. 

In the absence of the magnetic field we may set J = since the true ground state has 
zero angular momentum. In this limit, the above equation reduces to 



2 z 7 



9 E 



i 



(17) 



which implies that the confinement energy is v times the interaction energy. For v — 1 both 
are equal, which is not surprising, since in this limit the interaction energy scales like the 
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kinetic energy term. We may in fact regard the above statement as a "virial" theorem valid 
for all equilibrium configurations. 

We have, for the energy at the extrema, 

E = l i [(l + a 2 ) ( j ) + 2aJ + ^-]+g £ -4- 



2 1 
2u 



l + a 2 



aJ + 



v-U 2 
2i> ~6' 



(18) 



where <fi = R 2 <p and <fi is independent of J. This is the energy of the equilibrium configuration 
in a given J sector. A few comments are in order here. Consider the case when g = 0, which 
corresponds to the case of non-interacting confined particles in a magnetic field. This is an 
exactly solvable problem. ^From eq.([H]) we have, 



;i + a 2 -^)0 = O, 



(19) 



which implies either (1) = and therefore J = or (2) (1 + a 2 )<f) 2 = J 2 . The first of these 
solutions is the trivial solution which one obtains by directly extremising the Hamiltonian. 
The classical solution in this case is independent of the magnetic field. The second solution, 
on the other hand, can be obtained only by imposing the condition of fixed J during ex- 
tremisation. Equivalently, we vary the function F and not just H. The energy for a given 
J is then given by 



E 



(20) 



which is the energy of confined particles in a magnetic field. Further if we remove the 
confinement potential, the ground state energy is zero for all J < which is the solution for 
the lowest Landau level. (In the quantum mechanical case the zero point energy has to be 
added to the solution.). Thus the variational method in which the function F is extremised 
not only yields the correct energy but also the correct quantum degeneracy, which is infinite 
in this case. 

In the general case when g ^ 0, the second step in the extremisation involves extremising 
the energy with respect to J, that is dE/dJ = 0, 



[u +1)(1 + a 2 )R- 



va , J 



J 2 



9B? 



dR 

a/ 



2 R 2 



(21) 



since R depends on J. Differentiating R w.r.t. J in eq. (|T5|) and cancelling the overall powers 
of R, we have 



(z/ + l)(l + a 2 )i?-(z/-l) 



J 2 

PR 3 



dR 



J 



dJ 6m 2 ' 



(22) 
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Eliminating the derivative term in the above equations, we get the following solutions at 
equilibrium for each configuration: 



J = -a<pR 2 ] E = V -^~^- (23) 



Substituting for J in eq.fll5D, we get 

R 2 = [iguA(u,N))^, (24) 

where 



JV-l 



A[V - N) = g (1 + .J- 2., (25) 

An important point to note here is that this energy E is independent of the magnetic field 
and its dependence on g is explicit. The dependence on N and v is, however, involved. The 
angular momentum J extremising the energy depends on the magnetic field and is zero in 
the absence of the magnetic field as it should be. The expressions given thus far, though, are 
valid independent of the geometry of the clusters and are exact (for approximate solutions 
see eqs.(8,9) in ref. [§] for the special case of Coulomb interaction). 

We now specialise the general results given above to specific configurations which also 
happen to be ground states (global minima) for some N. While the results that we obtain 
are completely analytical and exact, the choice of configurations is based on the earlier 
numerical work ||. We have also checked the veracity of these results independently using 
numerical methods as outlined in Appendix B. 

The geometry of the clusters (or shells) is dependent on and A{y, N), which are as yet 
unspecified. The equations for the equilibrium configurations admit many solutions for a 
given N and v. In particular, there are two special configurations which are always solutions, 
viz (i) all the N particles are on a circle, O an d (ii) N — 1 particles are on the circle with 
one particle at the center, 0. For these two cases, only the overall scale factor R is to 
be determined. The angles %/2 are simply multiples of tt/N and ir/(N — 1) respectively. 
These configurations, however, need not be local minima in general. In particular, it has been 
numerically proved that for N < 5 the circle configuration is indeed the global minimum 
energy configuration (ground state) whereas for 6 < N < 8 it is the circle-dot which is the 
ground state in the case of Coulomb interaction ( v — 1/2 ). The ground state exhibits 
multiple shell formation for N > 9. In what follows we prove analytically that the first 
transition which occurs for N from 5 to 6 is independent of v. However, both O an d O 
remain ground states for iV = 5 and N = 6 respectively only for v up to a z/ min . In Fig.2, we 
display the eigenvalues of the Hessian (matrix of second derivatives) of the effective potential 
(see eq( |B.l| )) as a function of v. When iV = 5, the eigenvalues are positive definite up to 
v = 2.4 for both O an d O configurations. Both these configurations therefore correspond 
to local minima, but the calculation of the energy shows clearly that the O has lower energy 
and is in fact the ground state, confirming the earlier numerical calculations with v — 1/2. 
When iV = 6, the eigenvalues are positive all the way up to v = 10 (modulo a zero eigenvalue 
related to the rotational invariance) for the O configuration. Hence it is a local minimum. 
However the Q configuration ceases to be a local minimum for v > 0.484 and becomes a 
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saddle point since one of the eigenvalues becomes negative. It is possible that there are other 
configurations, other than the ones considered here, which are also local minima. They are 
not relevant to the analysis that is to follow. 





Figure 2: Eigenvalues of the Hessian of the effective potential as a function of v for 
O and O configurations and for N = 5, 6. Positivity of the eigenvalue 
indicates a local minimum. 

The cases of circle and circle-dot are particularly simple since there is only one scale 
involved, that is, all sf = l,i — 2,...,N and s\ — \ for the circle and s\ — for the 
circle-dot. 

For the circle case, we have, for iV particles, 

N-l 

<P=J2s 2 t + l = N (26) 



i=i 



and therefore the energy is given by 
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E Q = ^±A{4gvA N^]^ , (27) 



where 



N-l 1 1 N-l 



The second equality follows from the fact that for the symmetric configuration on the circle, 
9 Nk = 2n(N-k)/N. 

In the case of circle-dot, we have, for N particles , 

4> — N — l (29) 

since there are now iV — 1 particles on the circle and one at the centre. Therefore the energy 
is given by 

% = ^[^A Q (N - ir +1 ]^, (30) 

where 

A Q (u,N)=A (u,N-l) + l. (31) 

The extra 1 on the RHS is the contribution of the particle at the center. 

To ascertain which of the two configurations O an d O nas lower energy we look at the 
following ratio for the same number of particles: 

E Y +1 _ f n y +l a£> 



where 



N-l ^ 

fc=i sm (at) 

Note that in general the ratio / depends only on iV and v but not on the other parameters 
of the model Hamiltonian. Obviously the circle is a lower energy configuration iff / < 1. 
We claim that, for all v > 0, 

f(y,N) < 1 for N<5; 

f(u,N) > 1 for N > 6. 1 ' 

Further the function f(u, N) crosses unity exactly once for N between 5 and 6 and nowhere 
else. In Fig.3, we show the numerical values of the function f(u, N) as a function of ./V for 
various values of v. 
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N 

Figure 3: The geometric transition: f(v,N) vs N for z/ = 0.5 to v = 4 

The result can be easily checked for v = 1 since in this case A{J } = (iV 2 -l)/3. Therefore, 

N 2 (N+ 1) 



/(I, AO 



(AT — 1) (iV 2 — 2N + 24) ' 



(35) 



which reproduces the claims made above, for v — 1. For z/ << 1, again the proof is 
straightforward and we include it here since it will be relevant to the case of logarithmic 
interaction later. For very small v we make use of the identity that for any a, a u ~ 1 + 
z/log(a). Using this identity A^ may be written as, 



A 



N 



N-l 



(1 - 2z/log(sin(/br/A0) = N-l- 2uX N , 

k=i 



(36) 



where 



X 



N 



log 



n <) 



log 



N 



(37) 



Here we have used the identity, IlfcLi 1 sm (7f) — 2^1- Further we also approximate (N/N ■ 
l) u+1 « (TV/TV - 1)(1 + v\og(N/N - 1)). Substituting for \ N in /, we have 



N(N - 1 



(38) 



where 
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I* N = [N(N - 3) log(iV) -(N-1){N- 2) log(iV - 1)] 



(39) 



is independent of v. Clearly, whether / is less than or greater than unity depends on whether 
fiN is negative or positive. However the properties of /i/v cannot be dependent on v. It is 
now easy to see that /i^ is negative for N < 5 and positive otherwise. Hence the proof. 

Therefore, for N > 6, the ground state must have multi-shell structure (including possibly 
O). Note that the ratio / > 1 or / < 1 says nothing about whether the O or O is the 
ground state, or even local minimum. What it does say is that, if / > 1, then O can n °t be 
the ground state, independent of its stability properties. Therefore, whatever be the ground 
state it must have at least two shells (with one shell possibly trivial as in ©). This statement 
is independent of v. Thus we conclude that the first geometric transition is independent of 
v, g and J . Independent of the calculation of /, we know from numerical calculations that 
O is the ground state for 2 < N < 5 and O is the ground state for 6 < iV < 8 for v < u max 
where v max is determined from the eigenvalues of the Hessian. 

One could have also arrived at the conclusion that for N = 6 multi-shells must form, 
by studying the eigenvalues of the Hessian. This however is harder to compute analytically 
and also is not conclusive for small values of v for which O is still a local minimum. The 
criterion in terms of / > 1 is, however, easy and conclusive. 

It therefore appears that the organisation of many-body clusters in two dimensions into 
shells is a robust phenomenon, independent of the nature of the repulsive two-body interac- 
tion and also independent of the Hamiltonian parameters but dependent only on the number 
of particles in the cluster. In particular, the first geometric transition for the ground state 
from circle to circle-dot configuration occurs after iV = 5. The robustness of this transition 
seems to emerge purely from the number theoretic properties of the ratio of the energies in 
these two configurations. The results of this section are, however, restricted to parabolic 
confinement. In the next section, we analyse the effect of confinement on the geometry of 
clusters. 



III. EFFECT OF CONFINEMENT 



While the parabolic confinement is a simple and popular choice for model systems, real 
systems may be more complicated. In this section, we consider a generalisation of the form 
of the confinement potential and study its effects on the cluster properties. For our study we 



choose a form of the confinement potential proposed by Partoens and Peeters [16] who have 
studied the influence of the form of the confinement potential and the interaction potential 
using numerical methods as well as an approximate Thomson model. To keep the analysis 
simple we do not include the magnetic field. The Hamiltonian in dimensionless units may 
be written as 



N 



i=l L 



PS 
2 



~{r 2 V 
2 1 l> 



+ 9 £ 



1 



(40) 



where 7 characterises the confinement potential. While 7 = 1 corresponds to parabolic po- 
tential, we keep it arbitrary. The variational equations for extrema of F defined in equation 
(H) are straightforward to write. As before, one can eliminate the Lagrange multiplier in 
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favour of J and <fi (defined in eq.@). The equilibrium configurations are then determined 
by, 

J . , 

Pix = —rVi (41) 

Piy = -T X i ( 42 ) 

(^r'-g^^E^. (43, 

Taking the cross product in eq.(^3|) with fi gives N equations which are exactly the same 
as before in Section II. Taking the dot product with f and dividing by rf, however, gives 
the LHS which is dependent on the index i. Thus the parabolic confinement is special in 
this regard, since for 7 = 1 the LHS is independent of i. This fact was crucial for the 
conclusion that the shell structure does not dependent on the total angular momentum J. 
In the absence of the magnetic field, we may assume that the ground state has J = 0, in 
which case the z-dependence of the LHS may be easily eliminated as in Sec. I. The following 
analysis is therefore true for all J = equilibrium configurations in general and in particular 
for the ground state. 

Once again the shell structure is independent of the interaction strength. To see this we 
introduce an auxiliary variable, 

7 = EM) 7 - ( 44 ) 

i 

Expressing r*j = Rsi, R being a common scale factor (say the radius of the N th particle), we 
write 



(R 2 ) 



N-l 

2n 7 



E (O + 1 

L i=l 



i? 27 7 (45) 



since s% = 1 by choice, as in Sec. II. The basic equation determining the configuration space 
coordinates then takes the form 

Taking the scalar product with s, and dividing both sides by (s 2 ) 7 (7^ 0), we get 

, (R 2 Y +ry = y 1 - {Sj/Sj) cos(%) 

Note that the LHS is independent of the particle index i and the RHS is different from the 
parabolic case. However, we still have N — 1 independent of equations of the type 



E 



1 - (Sj/Sj) cos(%) 

)S- (7 ~ 1) (^ + ^-2^cos(^))^ 1 



1 - (sj/sk) cos(9 kj ) 
- k atT l) (4 + *? - 2a**i cos(^))^ 



E 2I n ^il^l^Ml , V M<- (48) 

1 2(7-1)/ 9 1 „2 "> ' 
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Further by taking the cross product with S{ and dividing by Sj, we get, 

V Sj S]n{9ij) = (49) 

which is independent of 7 and therefore is identical to the parabolic case. This equation 



provides a further set of N conditions on the internal coordinates S{. Together eqs. (fig 
provide the 2N — 1 necessary equations for determining the Sj and the angles Q{ and are 
completely independent of g and g'. 

These, the strengths of confinement and interaction potentials, determine the overall 
scale R through eq. (|47|) . Setting i = N in eq.([47f), we get 



, (R 2 y +J _ N yt l-3 jC os(0 Nj ) 
91 4gu fy(l + s]-2s jC os(9 Nj )y+i { ) 

Thus we have the result that the geometry or the shell structure of the equilibrium config- 
uration is independent of the parameters of the Hamiltonian, which only restrict the overall 
size of the system. 

Next we compute the energy of the equilibrium configurations. To this end we note that 
the auxiliary variable 7 defined in eq.(f44l) is related to the two-body potential energy by, 

g' 1 (j) 1 = 2gv ]T , (51) 

where the RHS is proportional to the potential energy due to interaction. Therefore, 

7[4£(n 2 ) 7 J = % E (52) 

which implies that the 7 times the confinement energy is v times the interaction energy. 
This is the generalised "virial theorem" valid for arbitrary confinement. We have, for the 
energy at the extrema, 

E = ^g'tr (53) 

This expression for energy is valid for all equilibrium configurations. The scale R may be 
calculated from eq. (|50D and is given by 

R 2 =[^A(u,N)}^ } (54) 
91 



where 

N-l 



^-g oT^i^^r - (55) 

and is independent of 7. This definition of A(u, N) is the same as in the parabolic case. 
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We now consider the effect of confinement on the geometric transition for iV from 5 to 6. 
As shown in the previous section for parabolic confinement, the configuration changes from 
circle to circle-dot independent of v. We first give expressions for the energy for these two 
configurations for any N. 

For the circle case, we have, for N particles, 



N-l 

[£4 7 + i] = iv 



(56) 



and therefore the energy is given by 



E 



O 



2v 



-9 



A Q N f 



g'l 



7 
i/+7 



(57) 



where Aq(u,N) is given in eq.fl28|). 

In the case of circle-dot, we have, for N particles, 

6^ = N-l 



(58) 



since there are now N — 1 particles on the circle and one at the centre. Therefore the energy 
is given by 



E, 



O 



z/ + 7 , 
-7S — 9 



91 w 



(59) 



where 



.4 



Q (u,N)=A Q (v,N-l) + 



1 



(60) 



as before. 

To ascertain which of these two configurations O an d O has lower energy we look at the 
ratio for the same number of particles, N, 



E 



o 



E, 



0, 



N 



A 



N 



N-l) + 

^(1-7) 



N 



N-l 



(61) 



where Ajy is given in eq. (|33D in the previous section. The effect of the confinement is 
therefore entirely contained in the prefactor whose exponent contains 7. 
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Figure 4: Effect of Confinement: /(7, z/, iV) vs A" for different values of '7. 
Tiie values of v range from 0.5 to 4 as in Fig. 3. 

In Fig. 4, we show the behaviour of / as a function of iV for typical values of v and 7. 
The pattern seen in Fig.4 may be summarised in two parts: 

1. For 7 < 1, the first geometric transition occurs before iV = 5. As v increases the 
transition moves towards iV = 5 but never beyond. In fact there exists v min after 
which the transition always occurs between 4 and 5 for N. 

2. For 7 > 1, the first geometric transition occurs after N = 6. As v increases the 
transition moves towards = 5. For v > u min the transition always occurs between 5 
and 6 for N. 



In the case of Coulomb interaction it was already known that the transition occurs either 
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before N = 5 (7 < 1) or after N = 6 (7 > 1) |L6 |. However the fact that the transition 
always occurs around N = 5 after some v min is, to our knowledge, a new result. 



IV. EFFECT OF THREE-BODY PERTURBATIONS 



In this section we consider a model Hamiltonian with a three-body interaction. While 
classically three-body interactions are not naturally introduced, quantum correlations can 
of course introduce effective many-body terms. The primary motivation for considering the 
three-body term here is to see if the shell structure survives under this perturbation. The 
effect of the three-body term is interesting and non-trivial on the first geometric transitions. 
The model Hamiltonian, in dimensionless units, is given by 



N 



i=i 



p, 2 M 1 2 

— r- 

2 2 1 



9i v ]_ 92 

+ 2 . f^,rf, + 2 



'ij ■ Tik 

<yi 2 ™,2 ' 
i,j,k(i^j^k) ij 



E 



(62) 



where g\ and g 2 are in general kept arbitrary. Notice that when g 2 = 0, this reduces to 
the model Hamiltonian analysed in section II with v = 1. The quantum mechanics with 
this model Hamiltonian (g 2 = 0) has some interesting applications to quantum dots and is 
considered in detail in ref. [TI|. The model Hamiltonian given above has a very interesting 
limit. When g\ = g 2 = g 2 , the quantum mechanical ground state and an infinite tower of 
states may be solved for exactly [|H],|1^]. Solving the Schroedinger equation, Hip Q = E ip , 
for the ground state we obtain 



^0 = n i^i 9ex p(-oE : 



(63) 



Kj 



and the ground state energy is given by 



E = N+-gN(N-l). 



(64) 



In what follows we keep g\ and g 2 arbitrary and comment on the case g\ = g 2 later. 

The classical analysis of the model Hamiltonian is done as before by extremising the 
Hamiltonian in the full phase space. As in the previous section we concentrate only on 
the case where the total angular momentum, J = , which implies that pi = and the 
equilibrium configurations are obtained as solutions of the equations 



251 E J-92 E 



r .ik 



+ 



Tjk 
2 „2 



.' jk ik ik ij . 



+ 292 E 

j,Kj+ k +i) 



Tjk-Tik -» Tik-^ij -» 
Tik ~r 7, 



2 4 
' jk' ik 



2 

r ik r ij 



■ (65) 



Again these equilibrium configurations may be local minima/maxima or saddle points. 

The general analysis of these equations proceeds as in section II. Since there are two 
coupling constants g\ and g 2 , the shell structure does depend on these parameters, unlike 
in the earlier models, except in the special case when the two are equal. Taking the dot 
product on both sides of the above equation with f\ and summing over the index i, we get, 



E^ 2 

i=i z 



9i 
2 



}_,92 



i,jU¥=i) % i 



r 2 2 ■ 
^ ijtjjtk 'ij'ik 



(66) 
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Therefore the confinement energy is again equal to the interaction energy including the 
two- and three-body terms. Using this, the energy in any equilibrium configuration may be 
calculated: 

N 

E = <p = R 2 Y J s 2 l =R 2 4>, (67) 
i=i 

where Sj are internal variables and R is the overall scale, which may be taken to be the 
distance of the farthest particle from the origin. We now specialise to specific configurations 
namely the circle and the circle-dot. It may be easily checked that these configurations are 
allowed equilibrium configurations for all N. These, however, need not be local minima for 
all N. For these two cases, only the overall scale factor R is to be determined. The angles 
Oij/2, as before, are simply multiples of n/N and ir/(N — 1) respectively. 
For the circle case, we have, for N particles, 

N-l 

0=[E^ + 1] = ^ (68) 

i=l 

and therefore the energy, after some algebra, is given by 

E Q = NR 2 = N^-^j^[ gi (N + 1) + 2g 2 (N - 2)] (69) 

In the case of circle-dot, we have, for iV particles, 

4> = N-1 (70) 

since there are now iV — 1 particles on the circle and one at the centre. Therefore the energy 
is given by 

% = (N - 1)R 2 = <N - U^p> + 2) + 92 (M^ + (N - 3),. ,7!) 

Interestingly, in the limit g± = g 2 = g 2 , we have the important result 

E = E = l -gN{N-l), (72) 

which is the same as the quantum mechanical ground state energy without the zero point 
fluctuation. The two configurations are also degenerate in this limit. We may therefore use 
either of these two classical configurations as a starting point for calculating the quantum 
corrections. Since the quadratic fluctuations reproduce the zero point energy the higher 
order corrections must be either vanishing or spurious. 

Coming back to the general case, g\ ^ g 2 , we would like to examine the effect of three- 
body terms in the interaction on the first geometric transition from circle to circle-dot 
configuration. Note that in the absence of the three-body terms, circle is the ground state 
up to 5 particles and circle dot is the ground state for N = 6. To ascertain which of these two 
configurations O an d O has lower energy in the presence of the three-body perturbations, 
it is sufficient to look at the ratio of the energies for the same number of particles: 
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= { EoY = (JL_) 2 N^[N + l + 2b(N-2)] 

1 Eq) \N-lJ (N- l)[N(N-2) + 24 + 2&(iV-3)(iV + 4)]' { ' 

where b is the ratio gijgx of the strengths of the two terms. Obviously the circle is a lower 
energy configuration iff / < 1. 
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Figure 5: f(N) vs N for various values of the ratio of the three-body to two-body 

strength, b = g 2 /g±. In (a) b varies from to 1 in steps of 0.1 and in 
(b) b varies from 1 to 2 in steps of 0.1. 



In Fig. 5, we show the numerical values of f(N) as a function of iV for various values of 
b = g^/gi- Some interesting facts emerge from Fig. 5. First, at b = 1, the ratio is 1 since 
the two configurations are degenerate in this limit. For b < 1, circle has lower energy for 
N < 5 and therefore / < 1. The ratio crosses unity when iV changes from 5 to 6 at exactly 
the same point no matter what the value of b is. It is as if the three-body perturbation 
has no effect on this geometric transition. The energies of course explicitly depend on the 
coupling constants. However, for b > 1, circle has higher energy for N < 5 and therefore 
/ > 1. Again the ratio crosses unity between 5 and 6 and for all N > 6, circle configuration 
has lower energy. This can be seen analytically also, from eq. (|73|) . 

However, it turns out neither the circle nor the circle-dot configuration can be the ground 
state when b > 1, that is when the three-body term dominates the two-body term in the 
Hamiltonian. In fact the ground state is one when all particles on a line, where the particle 
positions are determined by the zeros of the Hermite polynomial of order iV ||20|| . We 
must, however, point out that the above results are specific to the form of the Hamiltonian. 
Whether three-body interactions, in general, have a similar effect is a difficult question to 
answer. 
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V. DISCUSSION 



To summarise, we have discussed the shell structure of particle clusters in the presence 
of external confinement. To begin with, we have assumed harmonic confinement and the 
interaction is two-body. Following the analysis in section II, it appears that the organisation 
of many-body clusters in two dimensions into shells is a robust phenomenon, independent 
of the nature of the repulsive two-body interaction and also independent of the Hamiltonian 
parameters but dependent only on the number of particles in the cluster. In particular, the 
first geometric transition for the ground state from circle to circle-dot configuration occurs 
after N = 5. The robustness of this transition seems to emerge purely from the number 
theoretic properties of the ratios of the energies in these two configurations. A different 
type of confinement may in fact destroy the robustness of this transition for small v. The 
interesting point, however, is that after some critical u max , once again the results look similar 
to that of the parabolic confinement. 

Further, the presence of three-body terms in the interaction does not seem to affect the 
geometric transition as long as it remains a weak perturbation. Interestingly enough, within 
the model Hamiltonian we have analysed, there exists a critical point when the geometric 
transition changes dramatically as the perturbation strength is varied. Most importantly, 
the classical energy is the same as the quantum mechanical energy (without the zero point 
fluctuations) in the limit when the strengths of the two- and three-body terms are equal. It 
would be interesting to check if this is true of other many-body models where the quantum 
ground state is exactly known. 

While the existence of the shell structure is numerically well established, the analytical 
understanding of the existence of the shells is lacking. It would be interesting to know 
how many distinct extremal configurations there are . While the first transition seems to 
indicate the formation of shells, a deeper understanding of the Mendeleev table is clearly 
called for. It would also be interesting to study quantum corrections using these classical 
configurations as a starting point. Such a study would be of relevance to systems in which 
quantum corrections are non-negligible, as in, say, quantum dot systems. Some of the afore- 
mentioned questions are being probed further. 
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Appendix A: Logarithmic Interactions 

We briefly discuss the case of logarithmic interaction potential with harmonic confine- 
ment in this appendix. In some quasi-two-dimensional systems like quantum dots, this 
may be closer to the physical situation than the power-law interactions. Specifically the 
Hamiltonian, in scaled units, may be written as 



H 



N 



— = y 



-> 2 -Q. 

2 2 



9 E log hf , (A.l) 



where p is an arbitrary length parameter chosen such that p 2 > r 2 - V This is to ensure 
that the interaction potential is repulsive for all distances in the cluster. For simplicity, we 
do not include magnetic field here, though the analysis of Sec. I goes through identically in 
this case also. The conditions for the equilibrium are given by 

n = % E (A.2) 

Again, by setting = R 2 \£%^ s\ + 1] = R 2 (p, it is easy to demonstrate that the shell 
structure is independent of g (and also of the magnetic field, when included) as in Sec. 
I of the paper. Note that the equation above is manifestly independent of the arbitrary 
scale p introduced in the Hamiltonian. Further simplification occurs if we note that for any 
equilibrium configuration the auxiliary variable, is given by 

(j> = 2gN{N -1) (A.3) 

and the energy for the equilibrium configuration is given by 

E = gN(N-l)[l + \og(p 2 )]-g ]T log(rj). (A.4) 

We now calculate the energy for specific configurations, namely circle and the circle-dot. 
For the circle case, we have, for A" particles distributed symmetrically, 

E Q = gN(N - 1)[1 + log(p 2 ) - log(2£iV(iV - 1))] + gN(N - 3) log(iV). (A.5) 

In the case of circle-dot, for A" particles with N — 1 particles distributed symmetrically 
on the circle and one at the centre, the energy is given by 

E Q = gN(N - 1)[1 + log(p 2 ) - log(2£iV(jV - 1))] + g(N - 1)(JV - 2) log(iV - 1). (A.6) 

It is not possible in general to calculate these energies unless we specify the arbitrary 
scale p. However, to ascertain which of these two configurations O an d O has lower energy 
it is sufficient to look at the difference for the same number of particles, that is, 

AE = Eq — Eq = g[N(N - 3) log(iV) — (N — l)(N — 2) log(7V - 1)] = gp N . (A.7) 

where pn is already defined in Sec. II. Clearly the circle configuration has lower (higher) 
energy if p N is negative (positive). Note that this statement is independent of the arbitrary 
scale and also the interaction strength. It is now easy to see that p N is negative for A^ < 5 and 
positive otherwise. Therefore the first geometric transition from circle to circle-dot occurs 
when the number of particles changes from 5 to 6, as in the case of power-law potentials. 
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Appendix B: Numerical Simulations 

In this appendix we briefly discuss the numerical simulations carried out for multi-shell 
configurations. The analysis is carried out for the Hamiltonian given in Section II for 
parabolic confinement and may be adapted easily for other model systems discussed in 
Sections III and IV. Also collected are some results on the eigenvalues of the Hessian of 
the effective potential at the extrema. From eq.(^) it is easy to see that the same may be 
obtained as V14//(^) = 0, where 



V, 



eff 



;i + a 2 )0+^ 



+ 9 E j*. 



(B.l) 



The equilibrium configurations thus are obtained by minimising V e ff numerically. For this 
we use the version of the conjugate gradient method [plj which uses line minimisation. Since 
V e ff is positive definite and we use line minimisation, we are guaranteed to reach the extrema 
which are either local minima or some times saddle points but never local maxima. 

For a given number of particles, we choose several initial configurations graphically and 
compare the energies at the local minima. We also check that a local extremum reached 
is actually a local minimum by computing the eigenvalues of the Hessian of V e ff at the 
extremum. We verify the "Mendeleev" table of Bedanov and Peeters || completely for the 
special case of Coulomb interaction, i.e. when v = 1/2. 

The Hessian of V e ff can be computed easily. It turns out that of the 2N eigenvalues (for 
N particles), four can be obtained exactly. These also provide a check on the numerically 
determined eigenvalues. That the numerically determined eigenvalues contain the exact 
eigenvalues is also verified. Since the geometry of the extrema is independent of J, this 
analysis is done for zero magnetic field and J = 0. At an arbitrary point in the configuration 
the elements of the Hessian matrix of V e ff are given by 



d 2 V f 



eff 



dx;dxi 



d 2 V f 



eff 



dyidyj 



(1 + a 2 )5 ij -Agv 
(1 + a 2 )c% - Agv 



S H E 
S ij E 



y 2 k -{l + 2v)x 2 k . 



2 W+2 



x 2 k - (1 + 2v)y 



r s24-(l+2f)z 

Oij) 



r 



2 \u+2 



ik • 



r 



2 \v+2 
ik) 



;i + 2% 



r 



2 W2 
ijj 



(B.2) 
(B.3) 



d 2 V P 



eff 



d 2 V f 



eff 



dxidyj dyidxj 



k(k^i) V ik 



Vik^ik /-, r \ Vij-Eij 

2 \ v+ 2 _ V 1 _ ( r 2 y+2 



The eigenvalue equation for this matrix is 



9 2 Ve f! d 2 V eff 

dxidxj dxidyj 

d 2 Veff d*V eff 

dyidxj dyidyj 



Xi 



[l + a 2 )n 



X % 

Y; 



(B.4) 



(B.5) 



By grouping (Xi, Yi) as two-dimensional vectors Ri, one can write the matrix equation 
conveniently as 



1 + a 2 
Agu 



[l-n)Ri= E 



Ri 



2(1 + -) E V i = l...N. (B.6) 



ir 



2 \v+2 
ij) 
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At the extrema, r*j are of course solutions of the equations 

1 + a 2 \ _ ^ r, 
^ Tr 1 - 



i I n = E TTrrr V % = 1...N. (B.7) 



From these equations, four exact eigenvalues and eigenvectors can be deduced immediately: 

1. Both the equilibrium equations and the eigenvalue equations are manifestly rotation- 
ally covariant. One therefore expects p = to be an eigenvalue with the corresponding 
eigenvector Ri denoting the rotation of the r*j. This is indeed the case. Explicitly, p = 
and Ri = pkxfi V % — 1...N solves the eigenvalue equation when r*j satisfy the equilib- 
rium equations. Here p is an arbitrary, non-zero factor since the eigenvalue equations 
are homogeneous. 

2. If Ri = pri for all i = 1...N, then, using the equilibrium equations, it follows that the 
eigenvalue equations are satisfied with p = 2(1 + v). Thus the equilibrium configura- 
tions themselves constitute an eigenvector with an eigenvalue which depends only v. 
In particular it is common to all extrema for every N. 

— * 

3. It is immediately obvious that if R^ — a V i, then the RHS of the eigenvalue equations 
vanishes and p — 1 is the corresponding eigenvalue. As there are two independent 
choices for the two dimensional vector a, this eigenvalue is doubly degenerate. Note 
that this eigenvalue is totally independent even of the equilibrium equations. 

As all these eigenvalues are common to all extrema, the remaining eigenvalues must 
distinguish saddle points, local minima and local maxima. Apart from providing a check 
on numerical simulations, these should also prove useful in computing the semi-classical 
corrections to the classical energies. 
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